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We use density functional plus U methods to study the effects of a tensile or compressive substrate 
strain on the charge-ordered insulating phase of LuNiOa- The numerical results are analysed in 
terms of a Landau energy function, with octahedral rotational distortions of the perovskite structure 
included as a perturbation. Approximately 4% tensile or compressive strain leads to a first-order 
transition from an insulating structure with large amplitude breathing mode distortions of the 
NiOe octahedra to a metallic state in which breathing mode distortions are absent but Jahn-Teller 
distortions in which two Ni-0 bonds become long and the other four become short are present. 
Compressive strain produces uniform Jahn-Teller order with the long axis aligned perpendicular to 
the substrate plane while tensile strain produces a staggered Jahn-Teller order in which the long 
bond lies in the plane and alternates between two nearly orthogonal in-plane directions forming a 
checkerboard pattern. In the absence of the breathing mode distortions and octahedral rotations, 
the tensile strain-induced transition to the staggered Jahn-Teller state would be of second order. 


PACS numbers: 68.35.Rh, 71.70.Ej, 72.80.Ga, 73.61.-r 

I. INTRODUCTION 

The rare earth nickelates have been studied for many 
years m and have been of substantial recent interest 
IHai] following the proposal of Chalonpka and Khallinlin 
[ 6 ] that in an appropriately chosen superlattice configu¬ 
ration, an electronic structure similar to that found in 
the high-Tc cuprates could be realized. A one-band state 
has not been achieved, but the question of the degree to 
which the electronic structure can be controlled by ap¬ 
propriate combinations of strain and heterostructuring 
remains an area of active research lainHiaiiHiEQ]. 

The chemical formula of the rare earth nickelates is 
RNiOs, with R standing for La or for an element Nd, Pr, 
Sm, Gd, Eu, Lu of the rare earth series. The materials 
crystallize in variants of the ABO 3 perovskite structure 
with the R ion on the A site and the Ni ion on the B 
site. The basic structural motif is a corner-shared BOq 
octahedron. In the ideal perovskite structure the octahe¬ 
dron has six equal Ni-0 bond lengths and the point group 
symmetry of the Ni site is Oy. The important orbitals 
are the Ni-centered e^-symmetry d orbitals. Standard 
valence-counting arguments suggest that the Ni is in the 
low-spin dJ configuration with a filled t 2 g shell and one 
electron in the two e^-symmetry orbitals, which are de¬ 
generate in Oh symmetry. 

Having a single electron occupy two degenerate orbitals 
is expected to favor a symmetry breaking distortion in 
which one of the Cg orbitals becomes preferentially occu¬ 
pied and the point symmetry of the Ni is lowered from 
Oh to D 4 /J,. However, such a distortion has not been 
observed to date in the nickelate materials. With the 
exception of the LaNiOs (which remains undistorted to 
lowest temperatures), the materials exhibit at low tem¬ 
peratures an ordered phase [T] characterized by two dis¬ 
tinct NiOe octahedra, one in which the six Ni-0 bonds 
are short (but approximately equal) and one in which the 
six Ni-0 bonds are long (but again approximately equal) 


[7]. This disproportionation is sometimes referred to as 
“charge ordering” ais] based on the idea that the ionic 
charge of the Ni with longer Ni-0 bond lengths should be 
larger than of the Ni ions with shorter Ni-0 bonds, and 
based also on a difference in size of measured magnetic 
moments between the two sites. Although the actual 
charge difference between the sites is very small miiiH], 
for simplicity we will refer to the disproportionated state 
as “charge ordered”. 

The charge ordering is at first sight surprising be¬ 
cause the dominant interaction in transition metal ox¬ 
ides is generally believed to be a large on-site repul¬ 
sion “U” that acts to disfavor charge ordering. Indeed 
U = — 2 E^ is defined as the energy cost to 

change the electronic configuration from A'-electrons on 
each transition metal ion to the disproportionated con¬ 
figuration in which half of the ions have N 1 electrons 
and the other half have N — 1. The behavior is now 
understood as a consequence of a relatively large 

electronegativity of Ni. This places the rare earth nick¬ 
elate materials in or close to the “negative charge trans¬ 
fer gap” regime so that the electronic configuration is 
much closer to d^L than to d^: one electron is trans¬ 
ferred from ligand (oxygen) state to Ni so the Ni has 
two electrons in the eg orbitals (in the high-spin con¬ 
figuration) and there is an average density of 1/3 hole 
per O ion. Density functional plus dynamical mean field 
calculations [13 m] have shown that in this situation 
the “charge disproportionation” can be understood as a 
consequence of a hybridization (bond-centered) density 
wave leading to a site-selective Mott insulating regime. 
The high-spin d^ configuation of the Ni ions disfavors 
Jahn-Teller distortions with unequal occupancy of the 
Cg orbitals in agreement with measurements, indicating 
that even in strained superlattices of metallic LaNiOs the 
difference in occupancy of the two eg orbitals is small 
[ 121122 ] • Subsequent model system studies confirmed the 
essential features of this understanding m and suggest 
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that a negative charge-transfer energy implies that the 
effective low-energy theory is a two-orbital Hubbard-like 
model with a small or possibly negative effective U but 
a non-negligible J [20] . 

While this physical picture provides a satisfying un¬ 
derstanding of the essential features of the observations 
in terms of specific physics of the nickelate materials, it 
is incomplete in some respects. First, both experimental 
and theoretical studies of orbital disproportionation have 
focussed on the metallic regimes of the nickelate phase 
diagram [in na EiES] and leaves open the question of 
strain effects on the physics of the charge-ordered state. 
Second, the d^L configuration has the same symmetry 
properties as the configuration (this point was em¬ 
phasized by Peil et al m), meaning that the qualitative 
arguments suggesting a Jahn-Teller distortion should still 
apply. If a locally symmetric volume non-preserving dis¬ 
proportionation of the NiOe octahedra may occur, one 
may ask why not also cubic-tetragonal disproportiona¬ 
tions? 

In this paper we address these questions via an ab ini¬ 
tio study of the response of the charge-ordered insulating 
ground state of LuNiOs to an applied biaxial strain. In 
our study we fully relax the lattice subject to a constraint 
on the in-plane lattice constant, thereby approximating 
the effects of substrate-imposed strain on an epitaxially 
grown film. We find that strains on the order of a few 
% (i.e. of a magnitude comparable to those applied by 
epitaxial growth on reasonable substrate) have the po¬ 
tential to destabilize the charge-ordered state. For com¬ 
pressive strain the result is a metallic state with a mod¬ 
est cubic-to-tetragonal distortion of the NiOe octahedra, 
which is moreover approximately the same for each oc¬ 
tahedron. Sufficient tensile strain, however, is found to 
lead to a replacement of charge order by a spatially al¬ 
ternating in-plane Jahn-Teller order. We interpret the 
calculational results using a Landau theory free energy 
analysis which provides insights into the orders of the 
transition (indicating in particular that the transition to 
staggered Jahn-Teller state is intrinsically of second or¬ 
der and becomes first order only by virtue of competition 
with the charge-ordered state. Our results thus provide 
a different perspective on the strain control of orbital 
properties in transition metal oxides and show that the 
d^L configuration may also be susceptible to Jahn-Teller 
order. 


The rest of this paper is organized as follows. Section 
[TT| defines the system we study and the structural 
distortions we analyse and presents a Landau theory 
which encapsulates our results. Section in presents 


the specifics of our ab initio calculations. Section IV 
presents our computational results and Section [V| is a 
summary and conclusion. 
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FIG. 1. “Charge-ordered” structure of LuNiOs at vanish¬ 
ing external strain calculated using density functional plus U 
methods as described in Sec. |III| NiOe octahedra are indi¬ 
cated as grey cubes; the darker cubes have mean Ni-0 bond 
length smaller by O.lOA than the lighter ones. The triad on 
the left defines the lattice vectors. The calculated lattice con¬ 
stants |ai| = 5.I2A, |a2| = 5.52A, jaal = 7.36A are in close 
agreement with experiment [26]. We define Cartesian x, y 
and z coordinates so that z is parallel to as, and ai and a 2 
point approximately along the diagonals of the xy plane. 


II. FORMALISM 
A. Structure and strain 

In this paper we study LuNiOs. This material has a 
high charge-ordering transition temperature and an insu¬ 
lating ground state with a large gap to charge excitations 
[26] . The ground-state structure obtained from density 
functional plus U calculations described in Sec. |III| is pre¬ 
sented in Fig.[^ The unit cell has four inequivalent NiOe 
octahedra; in the absence of charge ordering, the octahe¬ 
dra differ only by rotations; the charge ordering creates 
two classes of octahedra with different mean Ni-0 bond 
lengths. Fig. also shows the lattice constants. From 
these we define a Cartesian coordinate system with ^ 
axis parallel to as and x, y axes in the plane defined by 
ai and a 2 but rotated by 45°. 

In the ground state of the actual material, the Ni-Ni 
distance in the basal [xy) plane is 3.76A, and there is a 
slight rhombic distortion, so the Ni-Ni bond angles are 
86° and 94°. We wish to simulate the effects of placing 
LuNiOs on a substrate, which will typically have a square 
symmetry. We therefore neglect the rhombic distortion 
and consider square structures with |ai| = |a 2 | and 90° 
Ni-Ni bond angles in the xy plane. We define the xy- 
plane lattice constant |ai| = |a 2 | = a as the diagonal of 
the square. In the undistorted structure, the lattice con¬ 
stant is a = a"*" ~ 5.3A, at which the energy is minimum. 
We will be interested in the consequences of a uniform 
compression or expansion of the lattice in the xy plane 
with the ^ direction free to adjust. We define strain 5a 
as an imposed change in the x^-plane lattice constant a: 


5a = a — . 


( 1 ) 
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The key variables in response to Sa are the shapes, sizes 
and orientations of the NiOe octahedra, which are the 
important dynamical variables of the structure. In the 
rest of this section we build up a theoretical description 
of the relevant distortions, starting from the simple case 
of an isolated NiOe octahedron and adding complexity 
as needed. 


B. An isolated NiOe octahedron 


To define notation we begin by consideration of one 
isolated NiOe octahedron. The unstrained structure is 
perfectly cubic (point symmetry Oh) with 6 mutually 
perpendicular Ni-0 bonds, which we take to lie in the 
±y and ±z directions. All the 6 bonds have the 
same length /q ~ 2A. 

The distortions of interest here preserve the inversion 
symmetry about the Ni and the orthogonality of the Ni- 
O bonds, so that a D 2 h symmetry is preserved. The 
distortions may be expressed in terms of three modes, 
defined in terms of the changes Sly, Slz in the x, y, 
and ^ bond lengths as 


Qo = 
Qi = 
Qs = 


Six T Sly + Slz 

Vi ’ 

Six - Sly 

V2 ’ 

-Six - Sly + 2Slz 


(2a) 

(2b) 

(2c) 


Here Qo is the volume expansion mode, Qi the (volume¬ 
preserving) x^-plane square-to-rhombic distortion, and 
Qs the (volume-preserving) cubic-to-tetragonal Jahn- 
Teller distortion in z direction. The Qo mode is invariant 
under Oh- The Qi and Qs modes together form a two- 
dimensional irreducible representation. Therefore, the 
energy of the octahedron ^(Qo,Qi,Q 3 ) to cubic order 
will be of the form 


structure at each Ni site. We must now attach a mo¬ 
mentum label to each mode. In addition, because the 
octahedra are corner-shared, there are constraints on the 
allowed momenta for each distortion. The momenta of 
interest are F = (0,0,0), G = (7r,7r,7r), C = (7r,7r,0). 
Note that these momenta are defined in the unit cell of 
the ideal cubic structure with one octahedron per unit 
cell. Of primary interest in interpreting the numerical 


results are the modes 

Qo = ; two-sublattice “charge order”, (4a) 

Qi = Q?; in-plane staggered Jahn-Teller. (4b) 

In addition, it will be useful to consider 

Qo = Qo 5 volume change, (4c) 

Qs = Qs ; uniform Jahn-Teller, (4d) 

Qs = Qs^; two-sublattice Jahn-Teller, (4e) 


which describe the response to a uniform strain and its 
coupling to a two-sublattice charge order. 

The energy function E'(Qo, Qs, ^s) of the 5 
modes is in general very complicated. A group theo¬ 
retical analysis is given in Appendix [A| The variables Qo 
and Qs are controlled by the strain Sa, which induces 
a Qs distortion and, via Poisson-ratio considerations, a 
nonzero volume change Qo of opposite sign to Qs- Both 
Qo and Qs are coupled to the order parameters qo, qi 
and qs, and these couplings will drive the phase transi¬ 
tions of interest. Our numerical results to be presented 
in Sec.|TV|may be understood in terms of an energy func¬ 
tion E{qo, qi\Sa) involving qo and qi only, with the other 
variables Qo, Qs and qs determined by the strain Sa and 
the values of qo and qi. 

Next, we focus on qo, the two-sublattice charge order. 
Viewed as a function of qo only (i.e. when qi =0), our 
results indicate that the energy E has a first-order tran¬ 
sition structure 

E{qo) = A2oql + ^ 40^0 + ( 5 ) 


E = AQl + B{Q\ + Ql) 

+ C (^Qi — Q 3 H-5 (3) 

where A, B and C are constants. We omitted the cubic 
terms Qq and Qo(Qi + Qs), which are just products of 
lower-order Oh invariants, and highlight the cubic cou¬ 
pling QiQs in the third term with coefficient C. In the 
lattice system, this part will give rise to an important 
coupling between the distortion Qf and the staggered 
Jahn-Teller order Qp, which we will define later. 

C. A corner-shared NiOe array 

We next consider an infinite 3D crystal of NiOe oc¬ 
tahedra, still with the Oh symmetry in the unstrained 


with strain-dependent coefficients ^ 20,^60 > 0 and 
A 40 < 0 near the transition. Thus, E{qo) has three local 
minima, at go = 0 and go = Eq'^. The strain Sa turns 
out to affect the value of q^ only slightly; the main ef¬ 
fect is to control via Qo and Qs the energy difference 
AE = ^(0) — E{±q'^), which is plotted against the in¬ 
plane lattice constant a in Fig. 

We next consider the energy as a function of gi only 
(i.e. go = 0). This energy will be found to have a second- 
order transition structure 

E{qi) = Ao 2 ql + Ao4qf, (6) 

with (at zero strain) A 02 , A 04 > 0 . Applying a strain Sa 
leads to nonzero Qo and Qs (see Fig.j^a)). Both of these 
may couple linearly to qf (recall in particular the QsQi 
cubic invariant), so that we have 

^02 = - Ag^Sa, 


( 7 ) 
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indicating that at tensile strains 5a > A^q 2 
plane staggered Jahn-Teller order could be favored. 

Finally, there is a biquadratic coupling A 22 Q 0 Q 1 be¬ 
tween the charge and Jahn-Teller orders; the sign of 
^422 > 0 is such that the two orders compete with each 
other. The resultant energy is therefore 

^cubic(^0: Ql) = E{qo) + ^(^l) + ^22^0^1* (8) 


D. Including octahedral rotations 

With Eq. (§ in hand, we now consider the structure 
of the actual materials. This involves a GdFeOs-type 
rotational distortion with four inequivalent Ni ions (see 
Fig. [^. The Oe octahedron around a given Ni site is 
rotated; the rotations may be symbolically written as 
Py , meaning that starting from the ideal cubic per- 
ovskite structure there is a rotation by angle a about the 
2 ; axis, and by angle /3 about the x and y axes. The su¬ 
perscript means the a rotations in neighboring oc- 
tahedra about the rotational axis of a (the 2 : axis) are in 
the same direction, while the ” means the p rotations 
in neighboring octahedra about the rotational axis of P 
{x or y axis) are in opposite directions. The angles a and 
P are small enough (< 15° in LuNiOs) that we may ne¬ 
glect the non-Abelian aspect of rotations and treat them 
as commuting (additive) axial vectors, rather than non¬ 
commuting (multiplicative) second-rank tensors. 

The important feature of the octahedral rotations is a 
breaking of the qi ^ —qi symmetry while preserving the 
qo —go symmetry. Here the k-points G = (tt, tt, tt) 
of go and C = (tt, tt, 0) of gi are defined with respect 
to the undistorted structure, i.e. with one octahedron 
per unit cell. This is allowed because even though the 
distorted and rotated structure now has 4 translation- 
ally inequivalent Ni ions, the Landau energy function of 
the system is still invariant under any translation by a 
nearest-neighbor Ni-Ni distance. 

The octahedral rotations generate an energy term that 
is linearly proportional to gi, and is of order a/3‘^ 10“^ 

in radians. The derivation of this is in Appendix us¬ 
ing group theory again. Similarly, a term Qoqi or Qsgi 
becomes allowed in addition to the Qsql term that we 
previously discussed. Thus, the final energy function is 
given by 

E = F/cubic ^01^1 

= F^cubic “ (^01^ + ^01 (9) 

where and are by a factor ofa/3^ 10“^ smaller 
than the coefficients in Fcubic? and Aq^^ results from the 
coupling terms Qoqi and Qsgi- The linear term is found 
to change sign at a compressive strain 5a = — Aq^^/Aq^\ 
We will look into the details in Sec. nvi The added term 
—Aoigi has an effect similar to that of an external mag¬ 
netic field on a system near a ferromagnetic transition. 


Eq. provides a minimal model that explains our 
data in all important qualitative aspects. In reality, there 
can be higher order terms of go and gi in both Fcubic 
and the symmetry breaking terms, as well as nonlinear 
dependence of the coefficients A^^ on strain 5a. 


III. METHODS 

In Sec. [TTj we constructed a Landau energy function of 
the bond-length distortion modes go and gi in KNiOs and 
took into account the GdFeOs-type rotational distortions 
/3y of the NiOe octahedra in a perturbative way. 
In this section and Sec. [TVl we study the structural tran¬ 
sitions of LuNiOs numerically by doing ab initio density 
functional + U (DFT+I/) calculations including struc¬ 
tural relaxation. The results are then interpreted using 
the Landau energy function in Eq. 

Our calculations use the Vienna ab initio simulation 
package (VASP) [27]. The LDA+I/ algorithm we choose 
in VASP is the rotationally invariant LSDA+V that fol¬ 
lows Ref. [28|. The Hubbard U of the Ni 3d orbitals in 
LuNiOs can be obtained by various methods, e.g. con¬ 
strained LDA [29l |30] , self-consistent linear response [31] , 
constrained RPA [32]|33], etc. They all give values of U 
within U = (5 ± 1) eV. The Hand’s coupling J is esti¬ 
mated to be 0.5 ^ 1 eV. We finally chose 1/ = 5 eV and 
J = 1 eV, as they gave a structure in Fig. that was 
closest to the experimental results. Slight changes of U 
and J within their errors were tried and no qualitative 
difference was found. 

We did a spin-polarized calculation using the PAW- 
PBE pseudopotential provided by VASP. The k-point 
mesh we used was 6x6x6 and the cut-off energy of 
the plane-wave basis was set to 600 eV. The errors due 
to k-points and energy cut-off are estimated to be smaller 
than the errors due to U and J by comparing with results 
obtained using a coarser k-point mesh of 4 x 4 x 4 or a 
lower energy cut-off of 400 eV. 

The computational unit cell was chosen to contain four 
LuNiOs formula units. Defining the basal plane as the 
one in which strain is applied, we take two formula units 
in the basal plane and two displaced vertically. To mimic 
the effects of a substrate, the in-plane lattice constants 
1 ^ 1 1 = |a 2 | = a are fixed to pre-set and equal values (so 
any in-plane rhombic distortion is neglected), jasj and 
all of the intra-unit cell degrees of freedom are allowed to 
relax. We slightly modified the conjugate gradient code 
in VASP to do this. 

The minimum energy of the substrate-constrained sys¬ 
tem is obtained at a = a"*" ~ 5.3A. The structure ob¬ 
tained is almost identical to the free structure in Fig. 
except that |ai| and |a 2 | are made equal (the small rhom¬ 
bic distortion is suppressed). We then adjust the sub¬ 
strate lattice constant a, our control parameter, away 
from a* and see how the structure changes. 
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(a) 
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FIG. 2. Upper panel: strain dependence of spatially uniform 
volume-changing (Qo) and even-parity volume-preserving 
cubic-tetragonal (Qs) octahedral modes. Lower panel: strain 
dependence of staggered volume-changing (qo) and two dif¬ 
ferent even-parity volume-preserving cubic-tetragonal octahe¬ 
dral modes {qi and ^ 3 )- Solid lines: results obtained from 
energy minimization. Dashed lines: results obtained from 
metastable states obtained by forcing staggered charge order 
(^ 0 ) modes to zero. 


IV. RESULTS 


A. Structures and energy difference 


Fig-i presents our main computational results: the 
evolution with strain of the structural parameters de¬ 
fined in Sec. |II C[ The upper panel shows the spatially 
uniform component of the relevant distortions. We see 
that an applied strain, as expected, induces both a uni¬ 
form distortion of the NiOe octahedra (the Qs mode) 
and a volume change. The changes are approximately 
linear in the applied strain. The lower panel shows that 
in the absence of strain the ground state is charge ordered 
{Qo 7^ 0)- Modest strain does not change the amplitude 
of the charge order but does activate a modest amplitude 
of staggered Jahn-Teller order (gs) as expected from the 
combination of strain and breaking of translational sym¬ 
metry. When the strain exceeds a critical value (in either 
the compressive or the tensile direction) the charge order 
vanishes via a first-order transition. On the compressive 
strain side, the resulting state is metallic (as seen from 



a/A 


FIG. 3. The energy difference AE = Ejt — Eco at differ¬ 
ent lattice constants a, with Ejt and Eco denoting the ener¬ 
gies of the metastable Jahn-Teller distorted structure (dashed 
lines in Fig.l^ and the stable charge-ordered structure (solid 
lines in Fig.^ between the transition points a ~ 5.1 A and 
a ^ 5.5A. Outside the transition points AE = 0 because 
the charge-ordered structure does not exist and relaxes to the 
only stable Jahn-Teller structure. 


the value of the density of states at the Fermi surface, 
not shown), and characterized by no order except that 
imposed by the strain. On the tensile strain side, the 
staggered charge order is replaced by a staggered Jahn- 
Teller order (gi) of comparable magnitude to the charge 
order. The Jahn-Teller order gi, unlike the charge order 
go, does not open a gap at the Fermi level, resulting in 
a metallic state. In this state, the long-bond direction 
of the Jahn-Teller order lies in plane. We also see that 
a very small-amplitude version of this order exists even 
in the charge-ordered phase, and the dashed lines show 
that if the charge order is suppressed, the amplitude of 
the Jahn-Teller order gi dramatically increases. 

The energy difference between the metastable and sta¬ 
ble states in Fig. is plotted in Fig. At zero strain 
a = a* ~ 5 . 3 A, the charge-ordered (CO) structure is 
lower in energy than the Jahn-Teller (JT) structure by 
82 meV per computational unit cell as defined in Fig. 
Under either a compressive strain (a < a"*") or a tensile 
strain (a > a*), the Jahn-Teller structure is favored and 
AE is reduced. At both transition points, the curve over¬ 
shoots a little bit to below zero, and ends at where the 
charge-ordered structure becomes locally unstable and 
relaxes to the Jahn-Teller structure. Both the overshoot 
and the linear AE — a relation near the transitions con¬ 
firm that the transitions are first order. 


B. The lower transition 

In this subsection we analyse the compressive strain- 
driven transition with the help of the Landau energy 
function in Eq. At fixed lattice constant a, we cal¬ 
culated the energy of a series of structures linearly in¬ 
terpolated between the structure with go = 0 and the 
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FIG. 4. The energy plots of linearly interpolated structures 
between the Jahn-Teller (qq = 0) and charge-ordered (mini¬ 
mum at qo = g*) states under compressive strains. The en¬ 
ergy of the Jahn-Teller structure with go = 0 is used as a 
reference point and the energies of other structures are mea¬ 
sured relative to it. The data points are fitted to Eq. §, 
with ^60 > 0 for all three curves. The other coefficients sat¬ 
isfy A 20 > 0,^40 < 0 for a = 5.115A and a = 5.125A, and 
A 20 < 0, A 40 > 0 for a = 5 .I 50 A. 


FIG. 5. Dependence of amplitude qi of staggered in-plane 
Jahn-Teller distortions on applied strain. Points are calcu¬ 
lated values. Solid line is result of fitting calculated points to 
Eq. (11). The solid line is the best-fit line and the dashed line 
is obtained by setting the linear coefficients = A^^ = 0 
in Eq. (11) to recover the ideal case of a second-order phase 
transition. The parameters of the best-fit line are = 

5.89 X 10“^, = 5.61 X 10■^ = 0.388, = 1.253, 

Ao 4 = 1, and a* = 5.30A. 


charge-ordered ground-state structure. Results are plot¬ 
ted in Fig.[^ We see that the energy has the typical first- 
order structure, with two locally stable minima crossing 
in energy as the lattice constant a is varied and we also 
see that the value of go characteristic of the charge order 
minimum is insensitive to the value of the strain. The 
main effect is simply a coupling of the strain to the en¬ 
ergy difference. 

The Landau energy function in Eq. © is reduced to 
Eq. ^ near the lower transition point a ^ 5.lA , because 
gi « 0 (see Fig.j^b)) makes go the only order parameter 
to consider. All coefficients A 20 , A 40 and Aqq are found 
to simultaneously change with the strain Sa. 

The transition can nevertheless be understood by a 
strain-induced change in A 20 alone. At zero strain a = 
^ 5 . 3 A (not plotted), the coefficient A 20 is negative, 
so that the ground-state structure is charge ordered. A 
compression of the in-plane lattice constant a favors the 
stability of the out-of-plane Jahn-Teller structure with 
go = 0, thereby increasing the quadratic coefficient ^ 20 - 
Sure enough, we see A 20 change sign in Fig. as a de¬ 
creases from 5 .I 50 A to 5 .I 25 A. Then go = 0 becomes 
lower in energy than qo = the charge order minimum, 
when a is further reduced. And finally, at some a below 
5.115A, the local minimum q = q^ disappears and go = 0 
becomes the only equilibrium structure. 

The transition from a nonzero qo = q'^ to 0 can be 
either first order or second order, depending on the sign 
of the quartic coefficient A 40 near the transition point 
a ^ 5.1 A. Since our data shows a first-order transition, 
we know A 40 < 0 and the charge-ordered phase with 
go = q'^ is stabilized by Aqo > 0 , i.e. the sixth-order 
term. This agrees with the fit parameters in Fig. 


C. Evolution of the Jahn-Teller structure 

The higher transition at a ~ 5.5A is more complicated 
because it involves both go, the breathing mode, and 
gi, the in-plane staggered Jahn-Teller mode, together 
with the GdFeOs-type octahedral rotations af f3y 
that break the gi —gi symmetry. We begin our analy¬ 
sis of this transition by considering calculations in which 
go is artificially set to zero. The heavy points in Fig. 
(equivalent to the dashed gi-line in Fig. [^b)) show the 
calculated evolution of gi with strain when go = 0. We 
see that over the whole range gi 7 ^ 0 , and that the evo¬ 
lution with strain is nonlinear. The nonzero gi is a con¬ 
sequence of the GdFeOs rotations, which, as previously 
discussed, couple linearly to the staggered component of 
the Jahn-Teller distortion. 

A minimal model to understand this evolution of the 
Jahn-Teller structure can be obtained by setting go = 0 
in Eqs. ©-(HI), leading to 

^(^1) = “ (Aq? + Ja)gi 

+ Ao2 - + ^ 04 ^ 1 , ( 10 ) 

where A 04 is assumed constant for simplicity. Eq. ( [1Q| ) is 
formally similar to the equation describing a ferromag- 
net in a magnetic field. The coefficients and A^^ 
are like an external magnetic field in the ferromagnetic 
case and arise from the breaking of gi — gi symmetry 
due to the GdFeOs rotations. The need to allow for a 
strain dependence of the coefficients is shown by the zero 
crossing of gi at a = ai = 5.20A. The dependence of A 02 
on strain reflects the tendency of tensile strain to favor 
the staggered Jahn-Teller order gi. 
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Minimizing Eq. (10) leads to 


+ 2 ^ 0 ? + 4 ^ 404^1 


^( 1 ) 

^01 




= Sa. 


( 11 ) 


We have fit Eq. (11) to the data points shown in Eig. 
and from the fit parameters we extracted the critical lat¬ 
tice constant a = a 2 = 5.61 (4)A at which the hypotheti¬ 
cal cubic structure would be unstable to staggered Jahn- 
Teller order, in the absence of charge order or GdEeOs 
rotations. We observe that while the uncertainties in¬ 
volved in fitting a four-parameter function to the data 
mean that individual coefficients cannot be determined 
with high accuracy, the estimated a 2 is robust. It is inter¬ 
esting that this value is not very much larger than 5.5A 
at which the charge order vanishes. 


D. The competition between qo and qi 


Comparison of the solid and dashed lines in Eig. 
shows that the staggered charge order (qo) strongly sup¬ 
presses the staggered Jahn-Teller order (gi). In the no¬ 
tation of Eq. the biquadratic term A 22 Q_qQ_i is large 
and repulsive. In terms of the analysis of Eq. (10) , ^02 
becomes ^402 + ^ 22^0 is so much more positive that 
until the charge order collapses at a first-order transition 
the staggered Jahn-Teller order cannot develop. There is 
therefore a strong competition between the two staggered 
orders go and gi. 


V. CONCLUSIONS 

We have used density functional and Landau theory 
methods to consider the effect of strain (induced by 
growth on a substrate with different lattice constants) 
on the charge-ordered state of LuNiOs. We find that the 
charge-ordered state plays a primary role in controlling 


the physics. It is the leading instability under ambient 
conditions, and its presence suppresses any other instabil¬ 
ities. However, with sufficient applied strain (within the 
DET+U approximation, of the order of 4%) the system 
undergoes a first order transition to a non-charge-ordered 
state. Interestingly, for tensile strain, the non-charge- 
ordered state is characterized by a staggered Jahn-Teller 
order. 

In the actual crystals, the symmetry breaking induced 
by the GdEeOs rotational distortion means that the stag¬ 
gered Jahn-Teller order does not break any additional 
symmetry of the system, but our Landau theory analysis 
indicates that even the ideal cubic nickelate will undergo 
a transition to staggered Jahn-Teller order if the tensile 
strain amplitude is sufficiently large. We thus conclude 
that even a system with a negative charge-transfer 
energy may have a Jahn-Teller instability. 

The actual magnitude of the strain needed to destabi¬ 
lize the charge order and allow other states is an impor¬ 
tant open question. While we imagine the strain as being 
produced by epitaxial growth on a substrate we have not 
included any quantum confinement effects in our mod¬ 
elling. Also, the DET+U method we have used is known 
to overestimate the tendency to charge order [21]. The 
charge order phase boundary also depends on how the 
double counting correction is implemented. More refined 
calculations, perhaps based on DET+DMET methods, 
should be employed to obtain better estimates for the 
strain needed to destabilize the charge order. But it is 
interesting that the magnitude of strain we have found 
is of the order of strains accessible by epitaxial growth 
on substrates. On the other hand, the first order na¬ 
ture of the transition means that there are no significant 
precursor effects to the transition. 
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APPENDIX 

A. The Landau energy function based 
on cubic (O/^) symmetry 


Once we find the symmetry-determined function E 
of the 9 modes, we project back to the 5 modes we 
previously started with by letting the other 4 modes 

= 0. The general form 

of E is then given by 


In the Fourier space, the nonzero distortion modes 
that appear in the calculated structures of LuNiOs (other 
rare-earth nickelates should have very similar perovskite 
structures) are Qo? Qs k-point F = ( 0 , 0 , 0 ), which 
are the uniform expansion and cubic-to-tetragonal dis¬ 
tortions of all octahedra, Qp at k-point C = ( 7 r, 7 r, 0 ), 
which characterizes a checkerboard pattern of in-plane 
staggered Jahn-Teller distortions, and at k-point 

G = (tt, tt, tt), which appear only in the “charge-ordered” 
structure. 

To find out the symmetry-determined form of the Lan¬ 
dau energy E’ as a function of the 5 modes 

and we need to extend our configu¬ 

ration space to a minimal Oh group-invariant subspace 
of 9 dimensions. This is because the Jahn-Teller dis¬ 
tortion along the 2 ; direction can be rotated to 

X and y directions by Oh to give us the mode. 

Similarly, rotating to x and y directions gives us 

And ~ can be rotated to 

gOTTTT ^ QnOn _ 

of course, mirrored to —and —We 
now see the 9 orthonormal modes 


/OOOO 
Vo 5 

er, 

nooo 
V 3 I 

/OTTTTTT 

Vo 5 

Qr", 


Qr°, 




(A.l) 


as a basis of the 9-dimensional extended configuration 
space that is invariant under Oh- 

The Land au energy as a function of the 9 modes 
in Eq. (A.l) will have to be invariant under the 3! = 6 
permutations of the x, and 2 : indices due to O/^, and the 
translations along x, y and 2 ; as well. A translation along 
X by one nearest neighbor Ni-Ni distance, for example, 
will leave all kx = 0 modes unchanged and let all kx = 
TT modes change sign. Translations in all 3 directions 
can generate totally 2^ = 8 ways of sign change. The 
Landau function E will therefore have to be invariant 
under 6 x 8 = 48 symmetry operations which include 
Oh + translations. At this point we forget about the 
boundary effects of the thin film and the substrate, so 
that the x, and 2 ; directions are all equivalent in the 
extended 9-dimensional configuration space. 

The algorithm we use for determining the symmetry- 
allowed form of the energy E is mainly based on the 
rearrangement theorem of group theory. We start with 
a general Taylor exp ansion of E with respect to the 9 
variables in Eq. (A.l) to some required order. The trun¬ 
cated expansion, which is a 9-variate polynomial, is then 
transformed by each of the 48 symmetry operations. The 
average of the 48 transformed polynomials is then guar¬ 
anteed to be invariant under all 48 symmetries according 
to the rearrangement theorem. 


00 2n 00 


^ = EEE C'njmiQO'! Q3 )Qo 


2n—j j 2m 

QsQi 


n=0 j=0 m=0 


(A. 2 ) 


where we have used the short-hand notations Qo = , 

Qs = QT. Qo = Qr", QS = Qi = Qr", which we 

also used in the main text. The functions CnjmiQo^Qs) 
are Taylor expandable and have the forms 

C'ooo(Qo5 Qs) = cio{Qo)Ql + bo{Qo, Q3)Qh (^-3) 

C'nOo(Qo, Q 3 ) = CiniQo) + ^n(Qo, Qs)^^: (^*4) 

^ 7210 (^ 05 ^ 3 ) “ ^ 72 (^ 05 ^ 3)^85 


where n = 1,2,3,''', and other Cnjm{QoyQ 3 ) func¬ 
tions and a ll lower-case functions that appear in Eqs. 
(A.3)-( A.5 ) are arbitrary Taylor-expandable functions. 
Eq. (A. 2 ) can be thought of as some advanced version of 
Eq. (3fin the main text with arbitrary constants aggre¬ 
gated into the Taylor coefficients of the arbitrary func¬ 
tions an{Qo), bn{Qo,Q 3 ), Cn{Qo,Q 3 ), etC. 

Now we study the strain effects, i.e. how things de¬ 
pend on the lattice constant a. Since the Qo and Q 3 
modes are at k-point F = (0, 0, 0), they are more closely 
related to the value of a than the other modes go 5 Q 3 i 
and gi. As a simplification, we assume that Qo = Qo(^i) 
and Qs = Q 3 (n) are smooth functions of the control pa¬ 
rameter a directly. As a increases, one expects Qo (<^)5 
the overall volume expansion mode, to monotonically in¬ 
crease and Q 3 (a), the overall cubic-to-tetragonal Jahn- 
Teller distortion, to monotonically decrease. The other 3 
modes go, g 3 , and gi may exhibit discontinuous jumps or 
other non-analytic behaviors at certain critical values of 
a and have to be kept as order parameters explicitly in 
the Landau energy function E. One may refer to the cal¬ 
culated structures in Fig. [^to see that the jumps in Qo 
and Qs at the transitions are much smaller. The energy 
function E is thus simplified to 


00 272 00 




72 J 772 




(A. 6 ) 


n=0 j=0 772=0 


which now has only three order parameters go 5 g 3 5 and gi. 
The modes Qo and Qs are treated as control parameters 
that are smooth functions of the lattice constant a and 
therefore disappear from the energy function. 

A further simplification can be made by noticing that 
in the calculated structures of LuNiOs (see Fig. [^, 
the order parameters go and gs, both at the k-point 
G = (7r,7r,7r), are always simultaneously nonzero, as in 
the “charge-ordered” structure, or simultaneously zero 
when the order is killed by a sufficiently large compres¬ 
sive/tensile strain. The fact that go and gs always go 
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“hand in hand” and “die together” suggests that we may 
combine them into one order parameter. This can be 
done by treating the ratio gs/g'o = ^ continuous 

(not necessarily monotonic) function of a. The Landau 
function is now further reduced to one with only two or¬ 
der parameters: 




(A.7) 


n=0 m=0 


where the coefficients 


2n 


-^ 2 n, 2 m (^) = E Cnjm{a)X^{a) (A. 8 ) 

j=0 

are arbit rary independent continuous functions of a. 
Eq. (A.7) gives the general form of the symmetry-based 


Landau energy function of RNiOs without considering 
perovskite octahedral rotations and non-orthogonal Ni-0 
bond angles. The Eq. in the main text is a simplified 
model that suffices to explain our numerical results. 


B. Including octahedral rotations: a 
perturbative approach 

Only even powers of qo and qi enter the Landau energy 
function E in Eq. (A.7). Here, we consider the effects of 


the perovskite octahedral rotations in RNiOs. We find 
that the rotations preserve the go ^ symmetry but 
break the gi —gi symmetry, which allows odd powers 
of gi to enter the energy function E. 

The rotational pattern in LuNiOs (and other rare- 
earth nickelates) is of the GdEeOs type, symbolically 
written as /3~Py . The meaning of the symbol is al¬ 
ready explained in detail in the main text. We want to 
construct an energy E as a function of the rotational an¬ 
gles /3“, Py and the bond-length modes go and gi 
using the symmetry group 04 /^ + translations, because at 
general values of the lattice constant a, the Jahn-Teller 
mode Q 3 (a) at the k-point F = (0,0,0) lowers the point 
group symmetry from (cubic) to (tetragonal). 

The benefits of using D 4 /J, instead of Oh are a) coupling 
terms involving Qo and Q 3 are automatically allowed for 
and b) the modes go and gi already form 

a group-invariant subspace without needing any exten¬ 
sions. Since all axial vectors /3“, and bond-length 
modes go, gi are invariant under spatial inversion I, only 
D 4 /i/{E, 1} = D 4 , which contains 8 symmetry operations, 
is effective in actually transforming the 5 modes. In addi¬ 
tion to D 4 , the translations can generate 4 possible ways 
of sign change according to the k-points of the 5 modes, 
among which and gi are at C = (tt, tt, 0 ), and Py 
and go are at G = ( 7 r, 7 r, 7 r). We therefore have totally 


8 X 4 = 32 symmetries to satisfy in order to construct a 
valid energy function E{az^ Px^ Py^ Qi) that takes into 
account the (small) octahedral rotations. 

Eollowing again the algorithm in Appendix [phased on 
the rearrangement theorem of group theory, we get the 
general form of the symmetry-allowed Taylor expansion 
of the energy function 


E —Aal -h B{PI -h Py) + Cgo + Dq\ 
T FdzPxPyQi T ' * * 5 


(B.l) 


where all coefficients A, H, C, D, F, etc. can be arbitrary 
functions of a. This is because the Qo(<^) and Q 3 (a) 
modes are functions of a and are invariant under D/^h 
and translations. They can therefore arbitrarily couple 
to any variables in Eq. (B.l). 


The omitted terms in Eq. (B.l) include other quar- 


tic terms that are products of the quadratic ones and 
higher-order terms. The leading-order term that breaks 
the gi —gi symmetry is EazPxPyQi = Eap‘^qi^ which 
is linear in gi. This term exists even if one considers the 
full Oh symmetry, which symmetrizes it to 


Fa^P^Pyqi F{azl3xPyQl'"° 

+ + OlypzPxQi^^)- 


(B.2) 


This means the octahedral rotations af P~Py break the 
gi —gi symmetry even at Q 3 (a) = 0 , i.e. at zero 
strain. However, the go ^ —go symmetry is strictly pre¬ 
served order by order. Switching the sizes of the larger 
and smaller NiOe octahedra in the “charge-ordered” 
structure is still a symmetry of the system even in the 
presence of the GdEeOo-type octahedral rotations. 

We therefore add the leading order symmetry-breaking 
term EazPxPyQi to the original Landau function E in 
Eq. (A.7) as a perturbation to get the symmetry right. 


The new Landau function is given by 


E = EE ^2n,2m{a)ql'^q1'^ + F{a)a^‘^qi. 

n=0 m=0 


(B.3) 


The added term should be small because aP^ ^ 1 for 
small rotations a and p. It should therefore be ineffective 
unless the even-power coefficients A 2 n 2 m{o) make gi = 0 
unstable or nearly unstable. The Eq. (j^ in the main text 
is a simplified model of the general Eq. ( |B.3[ ) here. 

Aside from octahedral rotations, non-orthogonal Ni- 
O bond angles can also break the gi — gi symmetry 
if the Ni-0 bond that is approximately along the 2 : di¬ 
rection forms different angles with the x and y bonds. 
The leading-order symmetry-breaking term should be 
also small and linear in gi, and can therefore be addressed 
in the same footing as octah edral rotations in our general 
model system given by Eq. (B.3). 








